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ABSTRACT 

We perform a linear analysis to investigate the dynamical response of a non- 
synchronized hot Jupiter to stellar irradiation. In this work, we consider the 
diurnal Fourier harmonic of the stellar irradiation acting at the top of a ra- 
diative layer of a hot Jupiter with no clouds and winds. In the absence of 
the Coriolis force, the diurnal thermal forcing can excite internal waves prop- 
agating into the planet's interior when the thermal forcing period is longer 
than the sound crossing time of the planet's surface. When the Coriolis effect 
is taken into consideration, the latitude-dependent stellar heating can excite 
weak internal waves (g modes) and/or strong baroclinic Rossby waves (buoy- 
ant r modes) depending on the asynchrony of the planet. When the planet 
spins faster than its orbital motion (i.e. retrograde thermal forcing), these 
waves carry negative angular momentum and are damped by radiative loss as 
they propagate downwards from the upper layer of the radiative zone. As a 
result, angular momentum is transferred from the lower layer of the radiative 
zone to the upper layer and generates a vertical shear. We estimate the re- 
sulting internal torques for different rotation periods based on the parameters 
of HD 209458b. 



1 INTRODUCTION 

Hot Jup iters are Jupiter-mass planets located within ~ 0.1 AU from their parent stars. 
Unlike Jupiter and Saturn in the Solar System, hot Jupiters are exposed to stellar irradiations 
that are much larger than their intrinsic fluxes. Consequently, a deep radiative outer layer 
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develops on the top of a convective interior in a hot Jupiter (e.g. see Guillot 2005 for a 
review) . 

Infrared observations of hot- Jupiter planetary systems with the Spitzer Space Telescope 
have been able to measu re temperature variations and therefore infe r temperature distr i- 



butions on hot Jupiters ([Harrington et al 



2006 



Knutson et al 



2007 



Cowan et al. 



20071 ). 



Meanwhile, a number of numerical simulations have been developed to investigate atmo- 
spheric circulation on a synchronized or non-synchronized hot Jupiter to better ascertain 
the origins of temperature distributions (see Showman et al. 2007 for a review). Despite the 
fact that these simulations are based on different equations and assumptions, and will thus 
exhibit different flow features, the simulated atmospheres usually end up with differential 
rotations such as banded structure or vertical shear. Although the flow patterns deviating 
from the initial uniform rotation are certainly the result of planetary rotation, the exact 
mechanism of how angular momentum is transported and redistributed between different 
regions of the atmosphere is yet to be established. 

When a global atmospheric flow follows non-synchronous rotation, the flow experiences a 
variation of stellar irradiation which serves as thermal forcing on the flow, producing thermal 
tides. Unlike ocean semi-diurnal tides which are driven by differential lunar gravity, the 



semi-diurnal oscillation of the atmospheric surface pre ssure 



on the Earth has been known 



to be mainly excited by the differential solar heating (IHaurwitz 



19641 ). In a state of quasi- 



hydrostatic equilibrium, the gravitational tide in the ocean and solid Earth and the thermal 
tide in the atmosphere can be modelled as gravitational and thermal bulges respectively (see 
Cartwright 2000 for a historical account). In the case of the Earth, the thermal bulge and 
the gravitational bulge have opposite phase difference with respect to the Sun (IHaurwitz 



1964 



Cartwrightl I200CH ) . meaning that the gravitational torques on the thermal bulge and 



on the gravitational bulge are pointing in opposite directions. Since Venus has a denser 
atmosphere and receives more solar insolation than the Earth, thermal tides on Venus are 
expected to be more prominent. This idea has inspired a number of models attempting 
to explain the slow retrograde spin of Venus by means of a balance between the torques 



due to gravitationa l and thermal tides (jGold fc Soter 



1969 



Correia et al. 



Dobrovolskis fc Ingersoll 



198(1 : 



20031 ). Laskar & Correia 2004 (cf. Showman & Guillot 2002) even postulated 



1 Diurnal tides of smaller amplitude also exist in the atmosphere at ground level (Chapman & Lindzen 1970, and references 
therein) but they correspond to a displacement of the centre of mass of the thermal bulge, which does not contribute to the 
gravitational torque (e.g., Correia et al. 2003). 
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that thermal tides may drive hot Jupiters away from synchronous rotation. This postulation 
suggests a mechanism of generating internal tidal heat in hot Jupiters and may lend support 
to the tidal inflation model (Bodenheimer et al. 2001; Mardling 2007 and references therein) 
in explaining why some of the transiting hot Jupiters are larger than indicated by current 
interior and evolutionary models. 

However, thermal bulges are probably not relevant to the case of gaseous (or liquid) 
planets. A perfectly rigid crust of a terrestrial planet can support any atmospheric pressure 
excess without being displaced sideways (or being slightly displaced if the crust is not per- 
fectly rigid; see Corriea & Laskar 2003). In the case of gaseous planets, the fluid underlying 
an overdense region is freely displaced sideways to attain hydrostatic equilibrium on the 
local sound crossing timescale. This means that any thermally driven density inhomogeneity 
on the top layer is almost cancelled out by the density inhomogeneity in the deeper layers^. 
By this argument, net thermal bulges cannot form on gaseous planets, and the gravitational 
torque acting on the thermal tide is essentially zero. 

Nevertheless, the oscillations of the stellar irradiation can still excite waves in gaseous 
planets. It is reminiscent of dynamical tides in the gravitational tide theories. Waves driven 
by gravit ational tides in hot Jupiters ha ve been studied in th e literature. Based o n the tidal 



Lubow et al. 



(119971 ) suggested 



theory by iGoldreich fc Nicholson! (119891 ) for high-mass stars, 
that the radiative layer of a hot Jupiter can be tidally synchronized by the internal waves 
excited resonantly by the tidal force of the host star. However, in contrast to high-mass 
stars where the external irradiation is unimportant compared to stellar intrinsic luminosity, 
the stellar irradiation onto a hot Jupiter is typically several orders of magnitude stronger 
than the intrinsic luminosity of the planet. It implies that the dynamics driven by stellar 
heating cannot be ignored. For instance, internal waves may also be excited thermally by 



stellar irradiation on the top of the radiative layer of a non- synchronized 



addition, rotation complicates the behaviour of internal waves. lOgilvie &: Lin 



r ot Ju piter. In 
(boOJ) studied 



the internal waves modified by Coriolis forces (i.e. Hough waves) in hot Jupiters. In the 
Earth's atmosphere, internal waves of the diurnal period are restricted in the region of 



2 One of the easiest ways to understand this concept is in terms of a planet covered by a liquid ocean and a gaseous atmosphere. 
If a thermal bulge is created in the atmosphere, then the surface of the ocean is displaced so that the column density perturbation 
in the atmosphere at each latitude and longitude is cancelled by an opposite column density perturbation in the ocean. In 
this way, the ocean can remain in hydrostatic equilibrium with no horizontal pressure gradients, because the same column lies 
above every latitude and longitude. 
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low latitudes where the Coriolis effect is small, and this explains why the thermal tide in 
surface air pressure i s 



1969 



Chapman Lindzenlll970l ). Semi-annual oscillations in Saturn's low-latitude stratospheric 



temperatures may be attributed to wave phenomena driven by seasonal thermal forcing 
(jOrton et al.ll2008l ). In the case of a hot Jupiter that is almost tidally locked by its parent 
star, the thermal forcing is much slower than the Coriolis effect, and will likely excite the 
Rossby waves (second kind of Hough waves; e.g. see Longuet-Higgins 1968). The importance 
of angular momentum transport by internal and Rossby waves has been discussed in the 
context of extrasolar giant planets (see, e.g., Cho 2008). It should be noted that while the 
waves driven by gravitational tidal forcing are able to exchange angular momentum between 
the planet and its host star, the waves driven by thermal forcing from the host star on the 
planet are only responsible for the angular momentum exchange between different parts of 
the planet, because of the cancellation of the gravitational torque described above. 



Atmospheric circulation is an extremely complex topic which involves turbulence, winds, 
as well as waves and how they are thermally driven and interact. Waves driven by thermal 
tides have never been studied analytically in the context of hot Jupiters to understand 
their basic behaviours. Therefore their roles in numerical simulations have not been easily 
identified. In this paper, we make a first attempt on the wave problem by considering a 
"clean" picture: a diurnal thermal forcing on the radiative layer with no clouds, winds, 
turbulence, and gravitational tides. The radiative flux in the atmosphere is modelled using 
the diffusion equation with a power-law Rosseland-mean opacity (cf. Dobbs-Dixon & Lin 
2008). Although the variation of the stellar irradiation is not small compared to its mean 
value, we employ a linear analysis and investigate the possibility of wave excitation in a 
non-synchronized surface layer of a hot Jupiter driven by stellar irradiation. The goal is 
to estimate how much angular momentum can be redistributed by thermal tides near the 
surface of a hot Jupiter in our simple linear theory. We first focus on the thermal tide 
problem for internal waves in a non-rotating plane-parallel atmosphere in §2. Then we turn 
our study to Hough waves in a rotating atmosphere in the form of a spherical shell in §3. 
Finally, the results are summarized and discussed in §4. 
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2 THE NON-ROTATING PLANE-PARALLEL ATMOSPHERE 

2.1 Basic equations 

We initially consider a non-rotating plane-parallel atmosphere with uniform gravity g = 
—ge z . The fluid equations for an ideal gas are 
du 1 

— + u-Vu = — Vp + g, 1 

at p 

^ + V-(pu) = 0, (2) 

dv 

-^ + u ' Vp + 7pV-u = -(7- 1)V-F, (3) 

where u is the fluid velocity, p is the gas pressure, p is the mass density, T is the tem- 
perature, F is the radiative flux density, k is the opacity, p is the mean molecular weight, 
R is the gas constant, a is the Stefan-Boltzmann constant, and 7 is the ratio of specific 
heats. For simplicity we assume that 7 and p are constant. We use the radiative diffusion 
approximation (j3J) throughout the atmosphere and apply the 'Marshak' boundary condition 
(cf. Pomraning 1973) 

aT 4 = \F Z + Fi (6) 

at z — +00, where F$ is the irradiating flux. The extension of the radiative diffusion approx- 
imation to the optically thin atmosphere is done for the sake of simplicity and is clearly a 
limitation of our model. 

2.2 Equilibrium state with a power-law opacity 

We consider an equilibrium reference state consisting of a static atmosphere that is uniformly 

irradiated by the mean stellar irradiation. For the equilibrium state we have 

dp 

16aT 3 dT 

F z = — = constant, (8) 

3np dz 

where F z is the intrinsic radiative flux density of the planet. Let r be the optical depth 
measured from z = +00. Then dr = —updz and we have 
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£ = (9) 

dr K 

^(^ 4 ) = |F,. (10) 
The solution of eq. (fTUjl subject to the boundary condition (jSJ) is 

ffT 4 = fF,(T + |)+F,. (11) 

Note that crT 4 = F 2 + F at the photosphere r = |. 

The equation for hydrostatic equilibrium can be analytically solved if we assume a power- 
law opacity: 

k = c K p a T- 4b (12) 



(13) 



for constants a, b, and c K . Then 
p a dp 4ag 
T Ab dT 4 = 3c^F' 
The solution satisfying p = at r = (where T = X^) is 

P a+1 = 4a 9 ( 1 A r T 4(fe+i) _ T 4(fe+i)l / u n 
a + 1 3c K F 2 U + l/ L °° J ' 1 J 

The top of the convective layer is located where the Schwarzschild criterion for marginal 

stability is satisfied; i.e., setting the Brunt-Vaisala frequency N = g 1 ' 2 [{l/^)d\a.p/dz — 

dlnp/dz] 1 / 2 equal to zero gives 
dlnp 7 
dlnT ~ 7-1 
at t = r conv . Thus 

1 conv\ _ JK 11 ' 1 ) _ v b+l 



(15) 



^ J 7 (a + l)-4(6 + l)( 7 -l) " { ' 

We require the denominator to be positive for convection to start. Since 

7^4 3 p f i 2\ , p 

1 _ 4 f ^' 2; 3^ Z J ' /i 7 \ 

T"4 if _|_ P. ' ^'^ 

-too 2 2 T 1 



we obtain 

4(F z /2 + F0 



(X-l). (18) 



If we treat k as a constant (a = 6 = 0), convection does not occur for 7 > 4/3. In this paper, 
the linear analysis will be performed for the radiative layer sandwiched by the top boundary 
at t = and the bottom boundary at r = r conv . 

Having found T(r) and p(r), we have p(r) and can then solve for z(t). However it is more 
convenient just to use r instead of z as a vertical coordinate in the problem. The solution is 
completely determined once the parameters g, c K , a, b, p, F z and Fj are specified. 
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2.3 Linear perturbation analysis 

We consider Eulerian perturbations of the form 
Re 



u'(z) e ik * x ~ iu)t 



(19) 



etc., where k x is a real horizontal wavenumber, x is the horizontal Cartesian coordinate, and 

u is a real frequency of the thermal forcing. In this paper, we shall consider a hot Jupiter in 

a circular orbit with the orbital period 27r/n or ;> and consider that its spin axis is normal to 

the orbital plane, although in this section we neglect the dynamical effects of rotation. The 

thermal tide is driven by a variation of the irradiating flux, and the problem at hand is to 

work out the amplitude and phase of the perturbations that result. 

The linearized equations read 
ik x p' 



-iuu x 

— iuju' z ■■ 

— \ujp - 

— iujp' - 
F' — F 

F' = F 7 



P 



~d z p' + ^d z p, 
P P 



u' z d z p + p(ik x u' x + d z u' z ) = 0, 

u' z d z p + ip{ik x u' x + d z u' z ) = -(7 - l)iik x F' x + d z F' z ), 
(ik x T'\ 



v d z T J 

' d z T' 3T' 
~cXT + ~T~ 



P_ 

P 



H 
K 



(20) 

(21) 
(22) 
(23) 
(24) 

(25) 
(26) 



l = Pi + r 

p p T 

This system of ODEs is of fourth order and the dependent variables can be taken as £ 2 , p', T' 
and F' z , where ^ z is the vertical displacement given by u' z = —iuj^ z . Rewriting d z = —npd T , 
we obtain the system 



d T £ z = (<9 T lnT - d T \np% 



k 2 x p' 



P 



UJ 2 Kp 2 



Kp \p 



T 



d T p -- 
d T T' 
d T F' 



k k \p T J ' 

v' 4T' F' 



d T T, 



P 



7 



,7 



d T In T — d T In p 



k 2 x F z T 



wp 

Kp 



P 
P 



7 



.7 



V 



(27) 
(28) 
(29) 
(30) 



K 2 p 2 d T T 

In the Appendix, we argue, using a scale analysis and a dimensional reduction of the problem, 
that the first term on the right hand side of eq. (128]) and the second term on the right hand 
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side of eq. (!30|) can be neglected. Neglecting these small terms amounts to assuming vertical 
hydrostatic balance and neglecting horizontal radiative diffusion. The large scales are also 
neglected since the geometry is planar and there is no rotation. 

The above four ODEs can be solved once four boundary conditions are given. In our 
model, we assume a thermal balance among perturbed energy fluxes at the top boundary; 
i.e., linearizing the Marshak boundary condition eq. ([6]) gives 

AaT 3 T' = -F' z + F[ (31) 

at r = 0. In other words, the thermal forcing, which is the perturbed irradiation is 
introduced to the system via the top boundary conditions. In the Appendix, we describe the 
mathematical details of how we determine the second boundary condition associated with 
the singular point at r = 0. 

To specify Ft, we assume that as the planet rotates, the stellar irradiation changes 
sinusoidally during the day and is completely switched off during the night. In the plane- 
parallel case, the stellar irradiation (heating term) is then proportional to 

cos H (cos 0). (32) 

where H is the Heaviside step function and is the longitude measured in a frame rotating 
with the orbit relative to the substellar point; namely, = x/R p — (n or b — Q)t. The above 
thermal variation can be decomposed into a Fourier series in as follows: 
11-2 

cos H (cos 0) = — I — cos H cos 20 + {m = 4 terms and above}, (33) 

7r 2 3n 

where m is the azimuthal wavenumber. The first term (i.e. m = 0) of the Fourier components 
is steady. It produces no tide but provides the uniform irradiating flux Fj. Other terms in 
the above equation give rise to the perturbed oscillatory irradiation F[. In this paper, we 
only consider the diurnal oscillatory component (i.e. m = 1, the second term on the right 
hand side of eq. fl33l ) for F[. In other words, the amplitude of F[ is {tz/2)F^ oj = n orb — Q, 
and k x — 1/R P for a planet of radius R p and spin rate Q. 

The other two boundary conditions at r = r conv are dependent on how the dynamics 
of the atmospheric gas varies with the thermal forcing frequency u (see the explanations 
following eq. (1371) for more details). Let c sp be the isothermal sound speed at the photosphere. 
Then the inverse of the horizontal sound-crossing time of the photosphere between the day 
and night sides of the planet is c sp /nR p . In the case of diurnal forcing, when \u\ ^> c sp /7cR p , 
the solutions of the ODEs behave like thermal diffusion and are expected to decay quickly 
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with depth. On the other hand, when \u\ <C c sp /ttR p , the equations admit a solution in the 
form of an internal wave which propagates downwards. If the depth of the radiative layer 
is large enough, the internal waves can decay quickly due to radiative loss before the wave 
reaches the turning point where N = Therefore in these two dynamical limits, we can 

set 



at r = r conv . We shall see later in this paper that the turning point of internal waves 
is extremely close to the bottom of the radiative zone. Hence when \u\ ~ c sp /ttR p and 
thereby enabling the internal waves to propagate to the turning point, these waves are not 
expected to have completely decayed at r conv . The boundary conditions p' — T' — may 
not be appropriate at r conv in this dynamical regime, and we should properly continue the 
wave solution into the convective region below. In this paper, we restrict ourselves primarily 
to the applications for large and small u> in the plane-parallel case. We note that setting 
p' — at the bottom boundary imposes the condition that the perturbed column density 
above the bottom boundary is zero as a result of the vertical hydrostatic equilibrium; i.e., 
this eliminates any thermal bulges in our calculations. As we have already explained in the 
Introduction, without a hard surface thermal bulges are unlikely to form on a gaseous planet. 

In summary, in the plane-parallel case we aim to solve the 4 linearized ODEs fl27j) - (!30j) . 
At t <C 1, we apply the boundary conditions associated with the imposed thermal forcing 
eq. (EH) (see eq. (1A18l) -eq. (1A21I) in the Appendix for details). At r = r conv , we adopt the 
boundary conditions eq. ( 1341) which are valid for \u\ ^> or <C c sp /ttR p . 

2.4 Numerical Results 

We consider an atmosphere with no clouds for simplicity in this paper, We choose HD 
209458b as an illustrative example for the thermal-tide study in the paper because the in- 
trinsic luminosity F z = 3.76 x 10 6 erg/cm 2 can be obtained from the simulation for the 
interior structure of HD 209458b in the grain-free case (Bodenheimer, private communi- 
cation). Some internal heating has been applied to the interior-structure model to explain 
the observed radius R p = 1.32Rj of HD 209458b (Bodenheimer et al. 2003). We also em- 
ploy the following input parameters for HD 209458b in solving the linearized equations: 
M p = 0.69Mj, Fi = 2.22 x 10 8 erg/cm 2 s, 2ir/n orb = 3.52474859 days, fi = 2 g/mol., and 
7 = 1.4 for the gas in the radiative surface layer of the planet. 



p' = T' = 



(34) 
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Figure 1. The structure of the uniformly irradiated equilibrium reference state for the radiative layer of HD 209458b. In the left 
panel, the data from Bodenheimer's simulation (cross points) based on the opacity table by Freedman et al. 2008 are compared 
to the analytical solution of eq. (1146 using the power-law opacity eq. (135 I t (solid curve). In the right panel, the relation between 
the pressure and the optical depth r derived from the analytical solution is shown. 



Without grains, we are able to fit the molecular Rosseland opacity computed by Freed- 
man et al. (2008) suitable for the radiative layer of a hot Jupiter (i.e. p ~ 1 bar and 
T « 1600K) to the power-law: 



where a = 0.42 and b = —0.987. c K = 10~ 17 ' 27 when expressed in CGS units. The fidelity of 
applying this power-law opacity can be justified by comparing with Bodenheimer's grain- 
free simulation. Figure [T] illustrates the comparison: the left panel shows that the structure 
of the radiative layer for HD 209458b from Bodenheimer's simulation agrees closely with 
the equilibrium state described by eq. (fT4"|) using the power-law opacity. We also apply this 
power-law opacity to the optically thin atmosphere r < 2/3. This is done for the sake of 
simplicity but is certainly a limitation of our model. We note that the radiative layer in the 
grain- free model is shallower than that in other interior-structure models (e.g., see Guillot 
2005). It is because R p and M p being the same, F z increases as the opacity in the radiative 
layer decreases, resulting in a thinner radiative layer according to eq. f[T8"|) . In the following, 
the vertical structure of solutions will be presented as a function of r. However, the readers 
can refer to the right panel of Figure [T]to convert the r coordinate to the pressure coordinate, 
for the case of HD 209458b. 

We set the positive x-direction as the direction of the planet's rotation. We focus on 
the case that the planet is rotating faster than its orbital motion. Therefore the thermal 
forcing and the thermal tides propagate in a retrograde sense; i.e., the substellar point 
moves backwards in the frame of the rotating planet (u < 0). 




(35) 
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Figure 2. Temperature perturbation T' as a function of r for the cases of 4 different forcing periods: —2iv/u) = 0.3, 1, 10, and 
100 days. The real and imaginary parts of T' are denoted by a solid and a dotted curve respectively. 
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Figure 3. Phase diagrams of T" in the complex-number plane for the cases of two forcing periods: 0.3 day (left panel) and 1 
day (right panel). The phase vector T 1 spirals clockwise from t sj to large t, meaning that T' exhibits a larger phase lag 
with respect to the star (i.e. the direction of the thermal forcing as shown pointing to the right in the horizontal direction) 
as t increases, its magnitude decreasing and finally dropping to zero. The phase vector T 1 for r ~ (i.e. at the top of the 
atmosphere) and r = 2/3 (i.e. at the photosphere) are shown. 



Figure [2] shows the vertical structure of T' in units of K for different thermal forcing 
periods ranging from 0.3 day (a case for a fast rotating planet) to 100 days (a case 

close to the synchronous state). The real and imaginary parts of T' are denoted by a solid and 
a dotted curve respectively. Although the results show that \T'\ <C T in all cases, \dT'/dz\ is 
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Figure 4. The Brunt- Vaisala frequency as a function of t resulting from the vertical structure with the use of the power-law 
opacity shown in Figure [T] In this case, convection occurs when r CO nv ~ 98.5. 

not smaller than \dT/dz\ especially at small r owing to the nonlinear forcing (i.e. \F[\ ~ Fi 
at the top boundary). Therefore our linear analysis is less justified. 

When the forcing periods are short (e.g. 0.3 and 1 day as shown in the top two panels), 
T' decays with depth and the solutions behave like those to the thermal diffusion problem 
with the heat diffusing from the top of the atmosphere to a depth characterized by the 
diffusion length rs ^2D/\u\, where D is the thermal diffusion coefficient of the atmosphere. 
Comparing the 0.3-day to the 1-day case, Figure [2] shows that T' can penetrate deeper 
in the 1-day result of a longer forcing period and therefore a longer diffusion 

length. The phenomenon of thermal diffusion can be also verified by the phase diagram of 
the complex number T' . Figure [3] shows that in the cases of the forcing periods = 0.3 and 
1 day, the solutions for T', denoted by the solid curve in the real complex plane, spiral 
clockwise toward the origin as r increases. The direction of the thermal forcing is shown 
by the horizontal arrow pointing to the right. For the forcing oc exp(— iuot) with u < 0, 
this means that the peak value of the perturbed temperature \T'\ exhibits a phase lag (i.e. 
delay) with respect to the star. Furthermore, while the phase lag increases with depth, \T'\ 
decreases with depth. All of these results demonstrate the process of thermal diffusion. 

On the other hand, when the forcing periods are long (e.g. 10 and 100 days as shown in 
the bottom panels of Figure [2]), the vertical profiles of T' exhibit wavelike solutions, meaning 
that waves are excited from the top of the atmosphere and propagate in. These waves are 
known as internal waves (i.e. g modes). The dispersion relation of g-mode oscillation in the 
WKB linear perturbation analysis without dissipation reads 
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2 2 N 2 2 

k x + k z = ~^2^x-> (36) 

where k x and k z are the horizontal and the vertical wave numbers respectively. Figure H] 
shows the vertical profile of N for our parameters for the radiative layer of HD 209458b. 
N ^> \u\ throughout the radiative layer except for the region very close to r = r conv « 98.5. 
The dispersion relation indicates that for a given k x , the vertical wavelength 2ir/k z decreases 
with the forcing period 2ir/u. Furthermore, internal waves can be dissipated due to radiative 
loss especially in the top layer of the radiative zone where the thermal timescale is short 
and k z is large. The shorter the wavelength is, the faster the wave is radiatively dissipated. 
This is exactly what is shown in the bottom panels of Figure O The internal wave for the 
100-day case has a much shorter wavelength and hence decays faster with depth than the 
wave for the 10-day case. 

The different dynamics appearing in the cases of short (thermal diffusion) and long (g 
mode) forcing periods can be also understood by comparing the forcing period with the 
sound crossing time of the planet's surface. When the forcing period is longer than the 
sound crossing time of the planet's surface; i.e., 

2. > , = ^ ^ Lg ,F±\ (*) ^ (tY 1/2 days , (37 ) 



M k x c sp c sp \1.32Rj J V1500K/ \2. 

downward-travelling internal waves (incompressible modes) can be driven by the diurnal 
forcing. On the other hand, when |u;| > k x c sp /2, the day-side and the night-side are causally 
disconnected^]. The diurnal forcing results only in the thermal diffusive effect in the vertical 
direction. 

In the plane-parallel case, internal waves carry a vertical flux of horizontal momentum. 
This flux builds up in the upper layers of the atmosphere where the waves are thermally 
forced, and is reduced at greater depth where the waves are radiatively damped and transfer 
their momentum to the atmosphere. It can be shown from the linearized equations (T2"Tl) - fl3"0|) 
that the gradient of the vertical momentum flux density is related to the non-adiabatic term 
V • F' as follows: 

rkr. 



rd T Re(pu' x u' z 



X -Re 



KpUJ 



h ^ (V • u')* (V • F') 



(38) 



lid 

Figure [5] shows the vertical profiles of the momentum flux Ke[pu' x *u' z ] and its gradient for 
the cases of two forcing periods: 10 and 100 days. The bottom two panels indicate that 

3 Although the day and the night sides are causally disconnected, the vertical hydrostatic balance of perturbations is still valid 
(see the Appendix). 
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X 



Figure 5. The top two panels show the momentum flux He[pu' x *u' z ] as a function of r for the cases of two forcing periods: 10 
(left panel) and 100 (right panel) days. In the bottom two panels, the gradient of the momentum flux multiplied by r, denoted 
by the dotted curve, is plotted for comparison with the solid curve described by the right-hand side of eq. J38P . 

the gradient of the momentum flux, which has been multiplied by r for clearer illustration, 
agrees with the right-hand side of eq. fl38l) . This validates the relation between momentum 
transfer and radiative damping in our result. General speaking, both cases indicate that the 
momentum flux is positive and indicate that the momentum is transported from the inner 
region where the gradient of the momentum flux is negative to the outer region where the 
gradient of the momentum flux is positive. In other words, the downward-travelling internal 
waves excited by the thermal forcing transport momentum outward. In the 10-day case, 
however, the situation is more complicated since the internal waves can penetrate deeper as 
a result of less damping due to longer vertical wavelength (see Figure [2]). 

We note that the undamped internal waves in the 10-day case reach the turning point 
where N = u (and therefore k z = 0) near the radiative-convective boundary (see eq. ( 1361) ). 
The waves are evanescent beyond the turning point and are reflected so as to interfere with 
the downward-travelling waves. As is clear from Figure HI the turning point is extremely close 
to the bottom of the radiative layer, implying that internal waves may not be evanescent 
significantly at r conv where we nevertheless have imposed the boundary conditions p' = T' = 
0. The consequence is that the solution in the 10-day calculation is sensitive to the bottom 
boundary conditions at T conv . To demonstrate this point, we consider the possibility that the 
perturbations at r conv may still preserve both the adiabatic and incompressible properties of 
internal waves instead of the "decaying" conditions p' = T' = 0. We apply alternative bottom 
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T 



Figure 6. Vertical profiles of T' (the solid and the dotted curve denote the real and imaginary part respectively) and the 
vertical momentum flux He[pu' x * u' z ] for 2tt/lu=10 and 100 days based on the following bottom conditions at T conv : p' /p = ■yp 1 / p 
and the Lagrangian density perturbation =0. 

boundary conditions to the 10- and 100-day cases: p'/p = ^p'/p (adiabatic perturbation) 
and the Lagrangian density perturbation =0 (incompressible perturbation). The results are 
shown in Figure [61 Comparing with Figure [2] and Figure [5], we find that the solutions in the 
100-day case are almost the same despite different bottom boundary conditions. However, 
the solutions in the 10-day case are indeed different for different bottom boundary conditions. 
Setting the bottom boundary conditions p' = T' = at a location well below r conv should 
give rise to more reasonable solutions for the 10-day case, but would require a method of 
treating the dynamics in the convective zone. 

The direction of momentum transport is related to the sign of u. A retrograde thermal 
forcing (u < 0) leads to upward momentum transport. When the sign of uj is switched to 
positive (i.e. the planets spins slowly than the orbit), our result shows that internal waves 
transport momentum downwards. 

3 THERMAL TIDES IN A ROTATING PLANETARY ATMOSPHERE 

We now consider the linearized dynamics of a thin spherical shell (a planetary atmosphere) 
which rotates at the uniform angular velocity fl We adopt the coordinates (8,<f),z), where 
9 and are the spherical polar angles and z is the altitude. In order to separate the vari- 
ables and determine the solutions, we consider the "perturbations" to refer to time- and 
azimuth-dependent deviations from a spherically symmetrically irradiated atmosphere. This 
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procedure corresponds to neglecting the latitudinal dependence of the average irradiation, 
and therefore eliminates winds in the basic state. Similar to the results for the plane-parallel 
case described in the preceding section, the damping of the vertically propagating waves 
should be able to transport and deposit angular momentum between different altitudes. 

3.1 Linearized Equations 

We adopt the linearized equations 

— iuou' d — 20, cos Qu' (j> = j^-dop', (39) 

-W. + 2Qcos04 = ^— r — j., (40) 

v pRpSmO 

= --d z p' + ^d z p, (41) 
P P 2 

- iup' + u' z d z p + pA = 0, (42) 



- iup' + u%p + 1P A = -( 1 - l)d z F' z , (43) 
d z T' 3T p' k'' 

H 

p p ~T 



imu' 



A = ^* ( *» sin «> + ^ + 9 ^ (46) 

where all perturbations have the form 

Re [u'(6, z) e im *-^] , (47) 

etc. We have assumed vertical hydrostatic balance and neglected horizontal radiative dif- 
fusion. These assumptions are identical to those made for the plane-parallel model and 
are justified in the Appendix. We have adopted the traditional approximation, which ne- 
glects 20 sin 9 u' z in the 0-momentum equation. The traditional approximation is valid if 
\u' z \ <C \u' d \, \u'A as expected in the atmosphere where the wave frequency is much smaller 
than N. Under these assumptions we can solve the horizontal components of the equation 
of motion for u' e and u'^ and substitute into the expression for A to obtain 

where £ is the Laplace tidal operator defined by 



Cp' = -tt(<9(9 + vm cot 9) 

sin 9 



sin 6* \ , 

(oe — vm cot 9)p 



I — v 2 cos 2 9 



sin 9 
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with v = 2Q/u and in general, u = m(n or b — Q). Together with regularity conditions at the 
north and south poles, this is a self-adjoint operator with real eigenvalues A (either positive or 
negative) depending on the dimensionless wave frequency v. In the non-rotating limit \v\ <C 1 
the eigenf unctions are associated Legendre polynomials (i.e. spherical harmonics) P™(cos#) 
and the eigenvalues are n(n + 1) for integers n > |m|; more generally the eigenfunctions are 
the Hough functions H v>m ^ n {6). It is traditional to express the eigenvalues in terms of an 
'equivalent depth' /iQ 

f-sH v,m,n X vmn H. vmri — H vmn . (50) 

We may then assume that all perturbations other than u' d and (which have been 
eliminated) depend on 9 through a particular Hough function. This allows for the separation 
of variables and we are left with the following system of ordinary differential equations: 

= --d z p' + ^d zP , (51) 
P P 

- iup + u' z d z p + pA = 0, (52) 

- icup' + u' z d zP + 7P A = -( 7 - l)d z F' z , (53) 
A = -^- + 5 X, (54) 

in which C has been replaced with the appropriate eigenvalue A. These are identical to the 
equations used for the non-rotating plane-parallel atmosphere and lead to the same ODEs 
(eq. (127]) - (130]) ) except that the horizontal wavenumber k x is replaced by the Hough eigenvalue 
A according to the formula k 2 x = X/Rp (or by the equivalent depth of the Hough function 
according to the formula k 2 x = uj 2 / gh) . In the non-rotating limit this means k 2 = n(n+l)/ R 2 , 
but in the rotating; be positive (wave solutions) or negative (evanescent solutions). 

A WKB analysis of the above linearized equations, in which (l/r)d o and d r are replaced 



by ike and ik r , gives the dispersion relation for adiabatic perturbations (jOgilvie fc Lin 



2004h 



V = K ( ^ ) , (57) 



uj 2 \ r 2 



where 



4 The reason for this name is that in Laplace's analysis of waves in a shallow incompressible ocean, the permissible values of 
v for free oscillations are determined by the condition that h equals the depth of the ocean. 
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A 



v 2 cos 2 9 



(58) 



This suggests that the solutions are oscillatory for A(l — v 2 cos 2 9) > 0. For solutions with 
A > 0, the oscillations are confi ned to the equatorial region. This correspon ds to the g-mode 



solutions modified by rotation (jLonguet-Higgins 



1968 



Bildsten et al. 



19961 ). However, if A is 



small, this confinement is relatively weak because ko is imaginary but small away from the 
equatorial region. The WKB analysis in the 9 direction starts to fail in this regime but the ra- 
dial (i.e. vertical) wavelength can remain small if w is very small. The solutions in this regime 



give 



rise to the baroclinic Rossby waves, or so called buoyant r modes (jLonguet-Higgins 



1968 : 



HeyJl200J). In the limit of solutions with very small, or zero, A, the radial WKB approach 



also fails. These special solutions that are global in both r and 9 are conventionally referred 
to as the barotropic Rossby waves, or simply the Rossby waves or r modes. Of course, our 
thin-layer calculation has eliminated barotropic r modes and we shall see later that non- 
adiabatic effects such as thermal diffusion affect the vertical structure of these modes to 
some extent. 

The above description of the solution properties implies that k x is not simply related 
to the horizontal scale of the tidal forcing, as in the non-rotating problem. There is a need 
for reinterpretation by finding the eigenvalues of the relevant Hough functions excited by 
the tidal forcing. The Hough functions for a given m and v can be deco mposed in terms of 



normalized associated Legendre polynomials P z m . (jOgilvie &: Linl I2004J ) show that the nth 



eigenvector of a certain tridiagonal matrix provides the nth Hough function in the form 

H v , m , n {9) = ]T a nJ Pr(cos(9)) (59) 
i 

for some coefficients a n ^ (the components of the eigenvector n). In the next subsection, we 
shall demonstrate how we determine the relevant Hough modes excited by thermal forcing 
and describe the associated problem with the method of separation of variables. 



3.2 Thermal Forcing 

Although the parent star is not distant from its hot Jupiter, we simplify matters by assuming 
the stellar irradiation to consist of parallel light rays impinging on the spherical planet. The 
stellar irradiation (heating term) is then proportional to eq. fl33l) multiplied by sin 9 with 
= 0- (n orb - Q)t. 

The latitudinal dependence sin 9 should be decomposed into Hough functions. If the 



Thermal tides 19 

Hough functions are expressed in a basis of associated Legendre polynomials with m = 1, it 
is necessary to decompose sin 9 similarly: 

sin^ = -(2/v / 3)P 1 1 (cos^) = -(2/^3) £ h >n H u>1;n (9) , (60) 

n 

where 6i >n is the first row of &z iTl which is the inverse matrix of a n> i in eq. (159|) . The coefficient 
&i jra tells us how much each Hough mode is excited by the latitude-dependent irradiation 
with m = 1 (i.e. the second term on the right hand side of eq. (|33|) ). 

Having said that the Hough modes are determined from the latitude-dependent heating 
in our model, we should note that our description is not entirely self-consistent because we 
have neglected the latitudinal dependence of the average irradiation. The ^-dependence of 
Fi would lead to a latitudinal variation of the properties of the unperturbed atmosphere and 
would also generate winds in the basic state. These complications are in conflict with the 
method of separation of variables used for solving the linear problem. However, as we have 
explained in the preceding subsection and as we shall see from some examples later in this 
paper, the Hough functions in some cases peak at low latitudes and decay quickly at high 
latitudes. Moreover, the latitudinal dependence of Fi, sin8, is not a fast varying function of 
low latitudes. These imply that if we aim to perform an order-of-magnitude estimate of some 
wave quantities integrated over all latitudes (such as the vertical angular momentum flux 
computed later in this paper), the contribution of calculations from high latitudes should 
be quite small. Therefore, we may simply apply the unperturbed heating at the equator to 
all latitudes and expect that the error introduced from high latitudes will be diminished by 
the Hough functions. This allows us to consider the "perturbations" referring to time- and 
azimuth-dependent deviations from the symmetrically irradiated atmosphere even though 
the excited Hough modes are still determined by the latitude-dependent perturbed heating 
in our model. 



3.3 Numerical Results 



We adopt the same input parameters for HD 209458b in the rotating case as in the plane- 
parallel case to solve the linear problem. We consider the diurnal thermal forcing (m = 1) and 
the scenario where the planet rotat es faster than i t s orbi t {u < 0). The Hough functions and 



eigenvalues are obtained based on lOgilvie &: Linl (120041 ) . We only consider the nth "wave" 
mode (i.e. A > 0)S which contributes the largest value of &i jTl ; namely, the largest heating 



5 However, A > does not necessarily admit a wave solution in the vertical direction in our problem involving thermal diffusion. 
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Figure 7. Hough functions H vm=ln=n i(0) = ^|Z^ 5 a n' ( cos f° r — 2-k/w = 0.3 (y ~ —2.17, solid curve), \ (y n —2.57, 
dashed curve), 3.5 iy fs —4, dotted curve), and 7 (i/ = —5.97, dash-dotted curve) days. 

term in the basis of Hough functions excited by the latitude-dependent diurnal forcing. This 
particular n is denoted as n! . Then we find that the positive eigenvalue A associated with 
the leading Hough mode for m — 1 has distinct features between fast and slow thermal tides 
(Longuet-Higgins 1968). 

When —2tt/uj is shorter than m 7.04 days (i.e. < 3 times the spin period or say —v < 6), 
A increases with the forcing period and is much larger than 1. The dominant Hough modes 
in this fast-tide regime have negative A and are evanescent. Positive and large A are then 
associated with less dominant Hough modes which normally consist of more weight from 
the associated Legendre polynomials of higher degrees. For instance, A for the forcing period 
= 0.3, 1, 3.5, and 7 days are 37.64, 53.9, 135, and 309 respectively. The corresponding eigen- 
vectors HvXn' i n the basis of the associated Legendre polynomials with the normalization 
coefficient \J[(2l + 1)(Z — l)!]/[2(7 + 1)!] summed from I = 1,3,5 up to 25 are depicted in 
Figure [71 These latitudinal structures agree roughly with the WKB analysis described by 
eqs. (JBTl) and (!58~]) : as the forcing period increases, the oscillatory solutions are more equato- 
rially confined and the latitudinal wavelength becomes shorter. These waves excited by the 
fast thermal tides are g modes modified by rotation. 

When the forcing period is exactly 3 times the spin period, A = and the solution of the 
Laplace tidal equations without thermal diffusion corresponds to a m — 1 Rossby wave with 
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Figure 8. Hough functions H u m ^ n=n i(ff) = J^iZi a n',l^l( coa ") f° r — 27r/w = 10 (v Ri —7.67, solid curve), 50 (v ~ —30.4, 
dashed curve), 100 (y fa —58.7, dotted curve), and 150 (f = —87.1, dash-dotted curve) days. 



the latitudinal profile oc sin 9(1 + 4cos 2 (0)). When the forcing period > 3 times the spin 
period, other than allowing solutions with a very large positive A (i.e. equatorially confined 
g modes), the tidal equations also admit solutions with a small yet positive A. They are 
known as the barotropic and baroclinic Rossby waves. These new solutions in the slow-tide 
regime are the predominant modes excited by the thermal forcing in our model. The A of the 
r modes increases slowly with the forcing period but remains smaller than 1. For instance, 
the A of the r modes for the forcing period = 10, 50, 100, and 150 days are approximately 
0.056, 0.109, 0.11, and 0.111 respectively. The Hough functions for these cases are plotted 
in Figure M to illustrate how the latitudinal structure of these dominant modes varies with 
the forcing period. Note that the r modes are less equatorially confined than the g modes 
and therefore couple better with the global heating profile (oc sin#). This explains why the 
r modes are more strongly excited than the g modes in our model. 

Knowing the eigenvalues A and assuming the unperturbed irradiation to be spherically 
symmetrical, we can solve for the z-dependence of the perturbations. For comparison, we 
start our study with the rotating cases for the same forcing periods considered in the non- 
rotating cases: 0.3, 1, 10, and 100 days. The results for the temperature perturbations 
are shown in Figure [9] (rotating cases) for comparison with Figure [2] (non-rotating plane- 
parallel cases). In the 0.3-day case, the solutions for both the rotating and non-rotating 
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Figure 9. Same as Figure [2] except that the Coriolis effect is included. 

cases exhibit diffusive behaviour. In the 100-day case, the waves propagate downwards in 
the rotating case as in the non-rotating case but the vertical wavelength in the rotating 
cases is longer (therefore the waves are less damped via radiative loss) because A has a small 
value. However in the 1-day and 10-day cases, the Coriolis effect introduces a large and a 
small A respectively, turning the diffusive solution to a wave solution in the 1-day case and 
turning the wave solution to a diffusive solution in the 10-day case. 

Although the leading modes in the 0.3- and 10-day cases give rise to vertically diffusive 
solutions, the less dominant modes (with therefore smaller magnitudes), do admit vertical 
wave solutions due to the much larger values of A. These are g modes and more equatorially 
confined. For instance, we find that (not plotted here) the second dominant mode in the 
10-day case gives a vertical wave solution (g mode) because A ~ 535 but its wave magni- 
tude, described by bi tU , is 0.024 which is much smaller than the magnitude bi }Jl t = 0.905 of 
the most dominant mode. The Hough function in the 10-day case has a latitudinal profile 
close to sin6 l (l + 4cos 2 (6 1 )), more akin to the barotropic Rossby mode. On the other hand, 
the dominant modes for the forcing period > 20 days do give vertical wave solutions. As 
the planet is more rotationally synchronized (i.e., —u becomes smaller and therefore —u 
becomes larger), the Rossby waves are more akin to baroclinic r modes with shorter verti- 
cal wavelength (see the 100-day case in Figure [9]) and are more equatorially confined. Note 
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that although the r mode in the 10-day case decays vertically due to thermal diffusion, it 
has a broad distribution across latitudes. Therefore the assumption of symmetrical unper- 
turbed irradiation relying on small values of Hough functions at high latitudes may not be 
appropriate in this case. 

In the regime of fast thermal tides {—v < 6), only g modes can exist. In the 1-day case, we 
have applied the bottom boundary conditions p' = T' = at r = 20 to illustrate the results 
(see Figure [H]) because we have difficulty obtaining the solutions when imposing the bottom 
boundary conditions at r = r conv . This may result from the possibility that the wave of long 
vertical wavelength in the 1-day case can penetrate deep into the radiative layer, rendering 
the bottom boundary conditions p' = T' = invalid (cf. the 10-day case in the non- rotating 
case). To justify this approach for the 1-day case, we examine the solutions for forcing periods 
slightly longer than 1 day, with the expectation that short-wave solutions would appear owing 
to the larger A. This is illustrated in Figure [TD] for the cases of —2tt/uj = 3.5 days (A ~ 135) 
and 7 days (A ~ 309). As the forcing period increases from 1 day to 3.5 days, and then 
to 7 days, the wave solutions can be seen clearly with decreasing vertical wavelength, as 
expected. The larger values of A associated with the less dominant Hough modes manifest 
themselves as shorter latitudinal wavelength as a result of the Coriolis effect, driving internal 
waves of shorter vertical wavelength according to the WKB dispersion relation (see eqs. ( 1571) 
and (EHD). 

In summary, the diurnal thermal forcing in the rotating case excites a series of Hough 
modes. When the forcing period is long (i.e., —2-k/uo is much longer than ~ 10 — 20 days), 
the dominant modes are baroclinic Rossby waves that can propagate downwards. When the 
forcing period is short (i.e., —2-k/uo is shorter than ~ 10 — 20 days), the dominant modes are 
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evanescent due to either the Coriolis effect (i.e. negative A) or thermal diffusion. However, the 
less dominant modes have short latitudinal wavelength and therefore internal waves can be 
excited without being subject to the same constraint that the sound crossing time between 
the day and night sides needs to be shorter than the forcing period as in the non-rotating 
plane-parallel case. Since the less dominant modes are small in magnitude, we expect that 
the angular momentum transported by the fast tides (g modes) is smaller than in the case 
of the slow tides (r modes). This is the subject of the next subsection. 



3.4 Vertical Angular Momentum Transport 

The vertical angula r momentum flux car ried by the Hough waves integrated over the planet's 
surface is given by (jOgilvie fc Lin! 120041 ) 

L = — f Re[p'*u' z }Rl sin 9d9, (61) 

LU JO 

where p' and u' z are both expressed as series of Hough functions (jp' = J2nPn^v,m,n-i e ^ c -) i* 1 
our model. Then the total flux is the sum of the contributions from each Hough function, 
because they are orthogonal. Therefore, the angular momentum flux carried by each Hough 
mode H v m n with m — 1 is given by 

L n = - T Re[p' n *u' zn ]H 2 u , m =i,nK sm 6dB = -Re[p' n * u' zn ]R 2 (62) 
to Jo to 

where we have used the normalization f£ H 2 mn sin 6d6 = 1. We focus only on L n due to the 
leading Hough mode (i.e., n = n') and use the peak value of the radial profile of the angular 
momentum flux to quantify the torque in each case of parameter study. 

Since our calculation is limited to a thin radiative layer and by the condition of symmet- 
rical unperturbed irradiation, to reasonably estimate L n >, we focus on the frequency regimes 
in which the solutions give short vertical wavelengths and modest equatorial confinement. 
Hence we carry out the estimates for the following forcing periods: 3.5, 7, 100, & 150 days. 
The results are listed in Table [U taking into consideration 4 different cases for each forcing 
period; Case I: the original case for the input parameters of HD 209458b, Case II: the same 
input parameters as Case I except that a larger k (c k is increased by a factor of 100) is used, 
Case III: the same as Case I except that a larger R p (= 2Rj, a young planet) is used, and 
Case IV: the same as Case I except that a smaller a (one half of 0.0474 AU) is used. Note 
that changing k, R p , or a would alter the interior structure and therefore other input param- 
eters such as F z need to change accordingly. However, the purpose of the present case study 



Table 1. Torques in units of 10 30 dyne cm generated by thermal tides with 2n/ui = 
cases. 
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-3.5, -100, and -150 days for 4 different 



Case 



—3.5 days —7 days —100 days —150 days 



I (original) 

II (high k) 



0.062 
0.022 
0.041 
0.056 



0.014 
0.005 
0.009 
0.013 



3474 
2276 
8414 
13719 



3011 
1973 
7294 
11893 



III (large R p ) 

IV (small a) 



here is simply to investigate how L n > varies with each parameter and to better understand 
what we can expect from our model for other interesting cases. 

In general, Table [1] shows that the torques in the fast-tide cases (i.e. 3.5 and 7 days) 
are weaker than those in the slow-tide cases (i.e. 100 and 150 days) by several orders of 
magnitude. As described in the previous subsection, this is primarily because the dominant 
Hough modes in the fast-tide regime are evanescent in the vertical direction. More specif- 
ically, the amplitudes of the internal waves |&i, ra '| are 0.08 for 2rr/u; = —3.5 days and 0.04 
for 2-k/lo = —7 days, which are smaller than the amplitudes of the baroclinic Rossby waves 
0.88 for 27r/c«; = -100 days and 0.82 for 27r/c«; = -150 days. 

Table [T] also illustrates how different input parameters affect the torque L n > driven by 
the thermal tides. Case II shows that increasing c K by 2 orders of magnitude only reduces 
the torque by less than a factor of 3, indicating that the result is less sensitive to the opacity. 
Case III shows that the torque driven by fast (i.e., 3.5 and 7 days) thermal tides on a larger 
hot Jupiter is weaker than that on a smaller hot Jupiter. The effect is opposite for slow 
(i.e., 100 and 150 days) thermal tides. Case IV shows that placing a hot Jupiter closer to its 
parent star decreases (increases) the torque in the fast (slow) tide case. Note that increasing 
n, R p (therefore g), or Fj reduces the pressure and density of the equilibrium state at a 
given t. This normally leads to an increase in u' z but does not dictate the change of p' in any 
unique way among the various cases. Therefore, the torque does not follow a simple trend 
in variation with any of the input parameters. 

Compared to the torque generated by the gravitational tides, the torque driven by the 
radiative damping of the baroclinic Rossby waves is sizable. With a constant tidal lag angle 
specified by the Q' value of the planet, the torque due to the gravitational equilibrium tides 
driven by a parent star on a non-synchronized planet is given by (e.g. Goldreich & Soter 



1966) 



L 



'grav 



3GM*R 5 p 
2a 6 Q> p 



(63) 
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which is comparable to the torques driven by slow Hough waves listed in Table HJ This 
implies that in our calculation, thermal tides play a significant role in transferring angular 
momentum to the atmospheric gas even if contributions of the gravitational tides are also a 
factor in the atmosphere. 

4 SUMMARY AND DISCUSSIONS 

We present a linear perturbation analysis for internal waves excited by the stellar diurnal 
thermal forcing (m — 1 thermal tides) in a non- synchronized radiative layer of a hot Jupiter. 
For computational simplicity, we employ the radiative diffusion approximation with a power- 
law opacity throughout our computation region from the radiative layer to a cloud-free 
atmosphere and apply the Marshak boundary condition for energy balance at z — +oo. We 
use the parameters of HD 209458b as an illustrative example. 

We first perform a linear perturbation analysis for a non-rotating plane-parallel atmo- 
sphere to explore how the dynamics of thermal response in the atmosphere varies with the 
thermal forcing frequency uj in the absence of the Coriolis effect. The boundary conditions 
at the bottom of the radiative layer are set to be p' = T' = which we speculate are less 
accurate for the waves with long vertical wavelength. Nevertheless, we find that when the 
thermal forcing period is shorter than the sound speed crossing time between the day and 
the night sides of the planet, the periodic thermal heating at the top of the atmosphere dif- 
fuses into a deeper layer with a decay length related to a;. In the fast thermal forcing regime 
corresponding to the rotation rate of our Jupiter, the problem for an atmosphere heated by 
the stellar irradiation becomes essentially a 1-D thermal diffusion problem with k x — > 0. On 
the other hand, when the thermal forcing period is longer than the sound crossing time of 
the planet's surface, the day and the night sides are causally connected and the problem 
becomes a 2-D phenomenon with incompressible properties. As a result, the internal waves 
are excited at the top of the atmosphere before propagating downwards. When the planet 
spins faster (slower) than its orbital motion, the thermal tides exhibit a retrograde (pro- 
grade) motion. The retrograde (prograde) waves causes the upward (downward) transport 
of momentum. The radiative damping of the waves leads to the deposition of the momentum 
in the atmosphere. 
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We then carry out a linear calculation for a thin spherical shell which rotates at a uni- 
form angular speed. Since |u;| <C N, we adopt the traditional approximation and neglect 
inertial terms in the vertical momentum equation. We also assume a spherically symmet- 
rically irradiated atmosphere as the basic state. These assumptions allow us to separate 
the variables and obtain the perturbation of the form 'Rje[u'(z)H vrrhn (0) exp (im(f) — iut)]. 
We examine the diurnal thermal forcing m = 1 and consider the nth component of the 
Hough function which contributes the dominant latitudinal structure of the heating. The 
linearized ODEs thus remain the same as those in the plane-parallel case except that the 
horizontal wavenumber k x is replaced by the eigenvalue y/X/R p . In other words, A conceals 
the information on how the latitudinal structure is modified by the Coriolis effect. Similar to 
the non-rotating plane-parallel case, the internal waves are driven by thermal forcing at the 
top of the rotating atmosphere. However, the internal waves in most of cases are primarily 
confined in a band of latitudes close to the equator and therefore are weakly excited by the 
global thermal forcing in our model. In the slow-tide regime (i.e., —2-k/uj is much longer 
than ~ 10 — 20 days), baroclinic Rossby waves can be largely excited by the global thermal 
forcing. When the planet spins faster than its orbit (u < 0), these waves propagate inwards 
but transport angular momentum outwards. While the torque generated by the radiative 
damping of the baroclinic Rossby waves in the slow-tide regime is comparable to the torque 
due to the gravitational tides, the torque generated by internal waves in the fast-tide regime 
(i.e., —2-k/uj is shorter than ~ 10 — 20 days) is smaller by several orders of magnitude. The 
magnitude of the torque is more sensitive to R p and a than k in our model. 



U nlike the thermal tide theories for a dense atmospher e on a terrestrial planet (IGold &: Soter 



1969 



Dobrovolskis fc Ingersolllll980l ; ICorreia et al.ll2003l ). our model for hot Jupiters, which 



do not have a hard surface, fails to generate net thermal bulges and therefore causes only 
internal transport of angular momentum inside the planet. In fact, the atmosphere in our 
model brings the interior closer to, not further away from, synchronous rotation. Together 
with the gravitational tides, one possible equilibrium state for rotation in our scenario is 
still the synchronous rotation. 

At this point, it is not clear how our linear results without winds for a non-synchronized 
hot Jupiter are able to explain the vertical shear appearing in existing 3-dimensional numer- 
ical simulations for the atmospheres of hot Jupiters. Apparently, the advection in the upper 
atmosphere of a hot Jupiter is driven by the temperature gradient between the day and night 
sides and cannot be described by our linear approach. Nevertheless, it is likely that the wave 
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dyna mics plays a role in th e lowe r atmosphere where the day-night temperature contrast is 



small. 



Cooper fc Showmanl (120051 ) postulate that the equatorial super-rotating jet occurring 



in the deep atmosphere of HD 209458b in their simulation could be due to wave transport. 
The numerical simulation by Showman et al. (2008) for a non-synchronously rotating HD 
209458b shows that the super-rotating jets are more equatorially confined and penetrate 
deeper in their case for 2k /oj = —3.5 days than those in their cases for 2tt/uo = 7 days 
and oo. Whether the deeper distribution of the super-rotating jet is due to the Hough waves 
driven by retrograde thermal forcing is an interesting subject worthy of further investigation. 

Our work presented in this paper is the first attempt at an analytical understanding 
of thermal tides in a hot Jupiter. Our linear theory suggests that the angular momentum 
transport in the atmosphere of a hot Jupiter due to a periodic thermal forcing is possible 
and mostly happens in low latitudes through the radiative damping of Hough waves. This 
encouraging result lays the foundation for further improvements in our calculations during 
future investigation, such as the inclusion of the semi-diurnal contribution, the considera- 
tion of Kelvin waves (third kind of Hough waves; e.g. see Longuet-Higgins 1968) driven by 
prograde thermal forcing, or the examination of internal waves driven from clouds in hot 
Jupiters. The extension of this work to other applications, for instance, thermal bulges on 
hot super-earths or Hough waves on a pseudo-synchronized hot Jupiter in an eccentric orbit 
(cf. Langton & Laughlin 2008), is necessary for gaining further insight into the question of 
how waves can play a role in atmospheric circulation on hot Jupiters. 
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APPENDIX A: DIMENSIONLESS EQUATIONS AND TOP BOUNDARY 
CONDITIONS 

In this section, we demonstrate how the fluid equations are non-dimensionalized in our 
problem. This procedure will not only help to clarify the relative importance among various 
terms in the equations more easily, but also shed light on the behaviour of the solutions near 
t = and help to specify the top boundary conditions. 

The equilibrium and perturbed states are completely determined once the parameters 
g, c K , a, b, /i, F z and F{ are specified. Therefore the problem can be non-dimensionalized 
in terms of these parameters. We can write p = pU p , T = TUt, p = pU p , k = hU K , and 
z — zU z , where p, T, p, k, and z are dimensionless functions of r and the dimensional units 
are given by the following relations: 

<tE# = F„ (Al) 
U P = yt , (A2) 

U p = -U P U T , (A3) 

u > = ctV' (A4) 

U K = c K U;U^ b . (A5) 



The equilibrium solution can then be expressed as 
f 



(^ + 1) + /1 1/4 , (A6) 



pa+l = _ 



4 /a + 1 



3 \b+ 1 



^4(6+1) _ ^4(6+1)1 ^ (A7) 



P = |, (A8) 
- = -/-^-' (A9) 

J pa T -Ab p v I 

where / = FJ F z is a dimensionless parameter that determines the importance of irradiation. 
A further dimensionless parameter associated with the problem is 

e= m (A10) 

where U c = (RUt/ p) 1 ^ 2 , a characteristic velocity unit. For the parameters adopted for HD 
209458b in this paper, e m 4.8 x 10^ 6 <^ 1. We further define an estimate of the radiative 
thermal diffusivity U x as follows 
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Uy = ■ (All) 

The linearized equations (!27I) - (l30l) can also be made dimensionless by writing lo = u}U u , 
K = k x U k , = l U z , p' = p' U p , T' = V U T and F' = P' F„ where 

U =^ 

w E/r 



E4 = ^. 

E/ c 

Note that UkU z = e. The above scaling leads to 

~ k 2 f)' 1 ( fi T' N 

= (d T lnT - <9 T lnp% - ^fW + — ' - 7^ 

Koo z p z up \p T 



d T p = -€ 
d T f'-- 

d T F' z =p 



K k \p T J ' 

v' AT' 
(a + l)?:-(6 + lW + ^ 
p T 



d T T, 



(A12) 
(A13) 

(A14) 
(A15) 
(A16) 



7 



,7-1. 



<9 r In T — <9 r In p 



• ~ % , 2 k 2 T 



k 2 p 2 d T T 



+ iup 

hp 



P 



7 



V 



(A17) 



v7- 1/ T. 

The terms proportional to e 2 may reasonably be omitted. Neglecting these small terms 
amounts to assuming vertical hydrostatic balance and neglecting horizontal radiative diffu- 
sion. 

For a + 1 > 0, we find the behaviour of physically acceptable solutions as r — > (omitting 
the e 2 terms) to be as follows: 

i = [C x + O(r)] lnr + [C 2 + 0(t)) , (A18) 
p> = t^tt [C 3 + 0(t)} lnr + r^r [C 4 + Q(t)} , (A19) 



T' 



C,r + 0(r 2 



lnr + [C 6 + 0(r)} 



F' z = t— [C 7 + 0(t)] In r + C 8 , 
together with the background states 
f = Too + Tir + O^ 2 ), 

p 



■•" +J = pn • 0(r 2 ), 



(A20) 
(A21) 

(A22) 
(A23) 



where = (| + J) 1 / 4 , 7) = ^(± + /)" 3 / 4 , and pi = (a + 1)(| + ff - The only variable 
to diverge as r is £ z , but since p = p/T = 0(r~), the mass flux at r = vanishes. 
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In fact the e 2 terms eventually become important as r — > 0, but we neglect them here. 
On substituting these series into the ODEs, we obtain the following relations between the 
coefficients: 



/ _ L2 T 2 \ | / 1 \"* 

C, = (. + !) ™ ^ C 3 (A24) 



Pi 



-46 

nh ~ a /( a + 1 ) \ T 



C. = -^-^s{*r\ ■ (A 26) 



■ooPi 

^ = ^#C 3 , (A27) 

Pi 

1 / 1 \" 4fc 

C 7 = -i^l /(a+1) d + (a + l)^— U- C 3 . (A28) 

Pi \Too / 

In addition, eq. fl3Tl) . the thermal forcing condition at r = 0, in dimensionless terms 
reads 

C 8 = 8T^C 6 - 2/', (A29) 

where the thermal forcing term /' = F(/F z . As a result, all the coefficients in the series 
expansion of the solution can ultimately be expressed in terms of C4 and C§. We apply 
the shooting method in solving the ODEs. We can first guess the values of C4 and Cg and 
evaluate the other expansion coefficients using the above relations. We then initialize the 
solution at r 1 and integrate to r = r conv , where two boundary conditions can be applied 
to determine the values of C 4 and C 6 and therefore the values of £, p', T", and F' at r -C 1. 



